Strength of effective Coulomb interactions in graphene and graphite 
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To obtain an effective many-body model of graphene and related materials from first principles 
we calculate the partially screened frequency dependent Coulomb interaction. In graphene, the 
effective on-site (Hubbard) interaction is Uoo = 9.3 eV in close vicinity to the critical value separating 
conducting graphene from an insulating phase emphasizing the importance of non-local Coulomb 
terms. The nearest-neighbor Coulomb interaction strength is computed to Uoi = 5.5 eV. In the long 
wavelength limit, we find the effective background dielectric constant of graphite to be e = 2.5 in 
very good agreement with experiment. 

PACS numbers: 73.22.-f, 73.22. Pr, 71.45.Gm 



The role of Coulomb interactions in graphene and re- 
lated materials poses a long standing problem: Experi- 
ments reported ferromagnetic ordering in nanographcne 
[H , in disordered graphite samples and at grain bound- 
aries in highly oriented pyrolytic graphite (HOPG) 
Ferromagnetism in pristine graphene, however, has been 
excluded experimentally for temperatures down to 2K 
Theoretically, the possibility of magnetism in de- 
fect free graphene has been predicted: An antiferro- 
magnetic insulating ground state has been obtained for 
the local Coulomb interactions exceeding a critical value 
Uaf > (4.5±0.5)t in Quantum Monte Carlo (QMC) cal- 
culations 0-0] and Uaf ~ 2.2i in Hartree-Fock theory 
[H, 0], where t « 2.8 eV is the nearest neighbor hopping 
parameter. A gapped spin-liquid has been predicted for 
on-site repulsion between U s \ = 3.5i and Uaf 0- Siz- 
able non-local Coulomb interactions can make the phase 
diagram even richer and lead to a competition between 
spin- and charge-density- wave phases 0, E3] or topolog- 
ically non-trivial phases IllJ. Doping of graphene might 
trigger further instabilities 1 1 2k 1 1 3| - In pristine graphene, 
the Coulomb interaction remains long ranged and it is 
controversial whether this might lead to strongly corre- 
lated electronic phases like an insulator 0, [3] or whether 
graphene is rather weakly correlated. The local part 
of Coulomb interaction is also crucial for the theory of 
defect- induced magnetism in graphene 15 1. 



The central issue in this discussion is the effective 
strength of the Coulomb interaction acting on the carbon 
jp z -electrons, which has only been estimated very roughly 
up to now [16[ : The bare on-site Coulomb interaction in 
benzene obtained from atomic carbon p z orbitals was es- 
timated to be 16.9 eV [13] • For polyacetylene, an analysis 
of optical modulation spectroscopy experiments within 
weak coupling perturbation theory yielded an effective 
on-site Coulomb repulsion of 10 eV [3, However, in 
this regime weak coupling perturbation theory might be 
inapplicable. For the long wavelength limit, reflectance 
measurements of graphite yielded a dielectric con- 
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TABLE I. On-site (Uoo, Uqo), nearest-neighbor (Uoi), next- 
nearest-neighbor (E/02, U02), and third-nearest-neighbor (f/03) 
(intra-layer) Coulomb interaction parameters for freestanding 
graphene and graphite. In graphene Uq% = Uq and Uq% = 
U 2 due to the sublattice symmetry. The bare and partially 
screened (cRPA) parameters are given. The cRPA parameters 
should be used in the effective Hamiltonian {JJ. 



stant of e = 2.4 due to screening by the high energy a- 

bands. This would correspond to an effective fine struc- 

2 

ture constant of a = w 0.9 for bulk graphite, where 
Tivf 



5.8eVA is the Fermi velocity 16|. For graphene, 
recent inelastic x-ray scattering experiments [2l| suggest 
a fully screened dielectric constant of e w 15 correspond- 
ing to a fine structure constant of a = 0.14. At the same 
time, first-principles GW calculations [22[ give e « 4, in 
agreement with the predictions of a simple Dirac model 
[16| . Recent experimental data on charge density depen- 
dence of the Fermi velocity [23| seem to be in agreement, 
rather, with the second value. So, up to now the strength 
of Coulomb interactions in graphene related materials 
has remained unclear and controversial — both theoret- 
ically and experimentally (for a review of correlation ef- 
fects in graphene, see Ref. |24I ). 

In this letter, we determine the Coulomb interac- 
tion strength in graphene and graphite within the con- 
strained random phase approximation (cRPA) [251 ] . We 
obtain ah initio effective Coulomb interaction param- 
eters that should be used in a generalized Hubbard 
model of graphene or graphite (see cRPA values in ta- 
ble UJ. We find that the on-site interactions in free 
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standing graphenc arc weaker than Uaf but close to 
the transition to the insulating spin liquid phase at 
U s i = 3.5t w 9.8 eV. Our calculations stress the im- 
portance of non-local Coulomb interactions in graphene. 
They put graphene in close proximity to two quantum 
phase transition lines and at the same time are possibly 
crucial for stabilizing a conducting state of freely sus- 
pended graphene. In the long wavelength limit, we find 
bulk graphite having an effective background dielectric 
constant e « 2.5, in agreement with the experiments from 
Ref. I20I . For graphene in the long-wavelength limit e is 
just one, as it should be for any two-dimensional system 
as will be explained below. 

We start with constructing a generalized Hubbard 
model for the graphene 7r-bands, 

H = -t ^2 c L c o^ - f c L c j^ 

<i,j>,o- <Cij»,rr 
i 

where o^ a annihilates an electron with spin a G {tii} 
at site i and n- h „ = c\ a Ci^. The index i = (i,A/B) la- 
bels the sublattice (A,B) and the unit cell centered at 
position Ri, f/y are the Coulomb interaction parameters. 
The nearest neighbor hopping is known to be t « 2.8 eV 
[Trj |26| and the next-to-nearest neighbor hopping t' de- 
pends on details of how the tight-binding parameters are 
determined: 0.02i < t' < 0.2t. 

To obtain all parameters entering the Hamiltonian (TTJ) 
from first principles, we performed density functional the- 
ory (DFT) and cRPA calculations. The DFT calcula- 
tions are carried out with the FLEUR code [27| using a 
generalized gradient approximation [28| for the exchange- 
correlation energy functional. We use a linear momen- 
tum cutoff of G max = 4.5 bohr -1 for the plane waves and 
an angular momentum cutoff of Z max = 6 in the muffin- 
tin spheres. The partially screened Coulomb matrix el- 
ements are calculated in the cRPA with the SPEX code 
29l 30| using the mixed product basis 29|, 31, 32| with 
cutoff values G' max — 4bohr _1 and L max = 4. 

The Hamiltonian (fT]) describes a system of C-p z elec- 
trons that interact via the effective interaction f/y, which 
incorporates the screening effects of all other electrons 
not contained in the Hamiltonian (J]). The cRPA ap- 
proach offers an efficient way to calculate this interaction 
[25j , as the screening channels are individually accessible. 
The two-dimensional symmetry of graphene clearly sep- 
arates the C-p z from other bands and, thus, enables an 
unequivocal elimination of the C-p z screening from the 
full RPA polarization function. Apart from the on-site 
term the resulting effective interaction yields the off-site, 
intra-orbital, and inter-orbital terms as well as their fre- 
quency dependence. 

The fully screened long wavelength dielectric constants 
reported in Rcfs. [2l| - |23| are different from the par- 
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FIG. 1. (Color online) Effect of lattice expansion on 
the strength of Coulomb interactions obtained with h = 
21.2 A. (a) On-site Uoo(a)/Uoo(a-o) and nearest neigh- 
bor Uoi(a) /Uoi(ao) Coulomb interaction as well as nearest- 
neighbor hopping t(a)/t(ao) as function of isotropic strain 
(1 — a /ao). a is the lattice constant. The parameters are 
given relative to their values at the equilibrium lattice con- 
stant ao = 2.47 A. The solid lines are linear fits serving as 
guide to the eye. (b) Ratios of on-site Uoo(a)/t(a) and nearest- 
neighbor Uoi(a)/t(a) Coulomb interaction to the nearest- 
neighbor hopping. The dashed line indicates the phase bound- 
ary at U B i/t = 3.5 separating the zero gap phase from a 
gapped spin liquid one @. 



tially screened cRPA dielectric constants obtained, here, 
in that the former include also contributions to screen- 
ing due to transitions between the graphene it bands. 
Hence, using the dielectric constants from Refs. 21-23^ 
in a generalized Hubbard model like Eq. (TTJ) or in the 
context of investigations like Refs. 0, [ljj would lead to 



double counting of screening terms arising from the tt 
electrons. 

We ensure the accuracy of the model parameters be- 
ing derived by carefully checking their dependence on the 
calculation procedure (the type of Wannier construction 
being used to define the C-p z orbitals) and convergence 
issues (Brillouin zone sampling and finite supcrcell height 
h) as we explain in the online supporting material 33 1. 
We find that Wannier functions directly from the C-p z 
projections [IH and 16 x 16 x 1 k-meshes for the BZ inte- 
gration yield accurate Coulomb interaction parameters. 

For graphene at its equilibrium lattice constant of ao = 
2.47 A, we obtain the Coulomb interaction parameters 
given in table Q] The on-site Coulomb repulsion Uqq B « 
3.3t is below Uaf « (4.5±0.5)t 0-Q but very close to the 
critical value of U s \ = 3.5i separating the zero gap phase 
from a gapped spin liquid one |8J. Comparing to the 
phase diagram reported in Ref. [lfj our results show that 
the nearest-neighbor Coulomb interaction of [/ 01 « 2.0t 
taken together with the local Coulomb interaction puts 
graphenc in close proximity to, both, a spin-density wave 
and a charge density wave transition line. 

The ratio of the kinetic energy given by t to the 
Coulomb interaction can, e.g., be changed by applying 
strain. Upon expanding the graphene lattice the near- 
est neighbor hopping decreases faster than the Coulomb 
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FIG. 2. (Color online) Static cRPA dielectric functions e(k) 
of graphene and graphite as function of (in-plane) momen- 
tum transfer k = (k x ,k y ). (a) Color coded (grayscale) e(k) 
for graphene. The most pronounced effect is the decrease of 
e(k) for k — ?> 0. There is a small directional modulation at 
intermediate momentum transfer k = |k| > lA . (b) cRPA 
dielectric functions e(k) of graphene and graphite as func- 
tion of k = |k|. Eq. ([3]) fits well the background dielectric 
screening for freestanding graphene in the limit of k — > 0. For 
graphite, two values of the perpendicular momentum transfer 
are considered: k z = and k z = 0.5(27r/c) with c = 3.3 A 
being the graphite interlayer spacing. 



interaction parameters (Fig. [TJa). An expansion of the 
lattice by a few percent leads to U~oo(a)/t(a) > 3.5, i.e. an 
increase of the ratio of local Coulomb interactions to the 
kinetic energy beyond the critical value of U s \/t = 3.5. 
In this situation, the non-local Coulomb interaction ef- 
fects can be crucial. It remains to be seen to which ex- 
tent the long range non-local Coulomb interaction screens 
the on-site repulsion (35| and stabilizes the semimetallic 
phase or whether non-local Coulomb terms drive the sys- 
tem towards other strongly correlated possibly topolog- 
ically non-trivial electronic phases as suggested in Refs. 

We now consider the Coulomb interaction in graphite 
and compare to graphene. In graphite, the two sublat- 
tices are not equivalent. We define the atoms of sub- 
lattice A be to directly above each other in adjacent 
layers and sublattice B as the atoms above hollow sites 
of the layer beneath. As table fl] shows, the on-site in- 
teraction in graphene and graphite is qualitatively sim- 
ilar with very little difference between the two graphite 
sublattices. The ratio of bare to cRPA nearest-neighbor 
Coulomb interaction is Uqi tc /U^ 1 ^ = 1.6 in graphene 
as compared to 2.2 in graphite. The non-local screening 
by the er-bands is considerably more effective in graphite 
than in graphene. 

This trend manifests clearly in the long wavelength 
limit as can be seen from the Coulomb interaction in re- 
ciprocal space. To this end, we consider the Coulomb 
interaction matrix elements in terms of the Bloch trans- 
formed C-p 2 -Wannier functions, |iy n k}- We calculate the 
ratio of bare to cRPA screened interaction [36j 
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FIG. 3. (Color online) Frequency dependence of the on- 
site and nearest-neighbor interaction obtained from cRPA 
for graphene (h = 21.2 A) and graphite. For graphite 
Uoo(u) — Uqq(uj) is shown, which is virtually the same as 
U&iJ). \u£(uj) - [/ s (a;)| < 0.15 eV for u < 20eV. 



For graphene, our cRPA calculations (Fig. [2]) yield 
e(k) ss 2.4 for intermediate momentum transfer, k = 
k| > lA -1 , and e(k) — > 1 for k — > 0. The screening 
due to high energy states in graphene becomes essen- 
tially negligible in the long wavelength limit |37j|. This 
is fundamentally different for graphite where e(fc||) *=s 2 
almost independently of the momentum transfer and 
e(k = 0) « 2.5. Hence, graphite should be less corre- 
lated than graphene. 

In the long wavelength limit, the simplest model to 
address screening by high-energy bands in freestanding 
graphene is to consider a film of thickness d and dielectric 
constant ei. Transferring Ref. [IH to the geometry at 
hand [3911 we obtain 



hH*) = - 



1 ei + 1 + (ei - l)e 



-kd 



ei + 1 



Ui 



kd[—-e 1 



h 1/2 



(3) 
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Our cRPA calculations confirm this expectation (see Fig. 
O. Eq. Q turns out to describe the partially screened 
Coulomb interaction well for /c = |k| < lA -1 with 
d = 2.8A and t\ = 2.4 proving the applicability of this 
classical model at long wavelengths. 

Integrating out the graphene cr-bands and other high 
energy states leads to frequency dependent effective 
Coulomb matrix elements. For graphene and graphite, 
the effective Coulomb interaction is significantly fre- 
quency dependent above u> > 5 eV (Fig. [3]) . Within the 
energy range of the Dirac spectrum, however, (~ 2 eV) 
the Coulomb interaction can be well considered in the 
static limit. 

In conclusion, the strength of Coulomb interactions 
in graphene and graphite is accurately determined by 
first-principles calculations. The local Coulomb inter- 

A/B 

action in graphene is J7 00 = 9.3 eV ~ 3.3£, which is 
very close to the critical value U s \ = 3.5i for the transi- 
tion to a gapped spin liquid. By straining graphene, the 
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system can be driven across this critical value. More- 
over, we find large non-local Coulomb interactions (e.g. 
Uoi = 5.5 cV ~ 2. Of). By means of a dielectric substrate 
below graphene the screening of the long range tails of 
the Coulomb interaction can be tuned, while the local 
Coulomb interaction terms are expected to be much less 
affected by the dielectric environment. Hence, also the 
ratio of local to non-local Coulomb interactions can be 
tuned. It remains to be seen which additional many body 
instabilities might be triggered in this way or to which 
extent the conducting state of free standing graphene can 
be stabilized by non-local Coulomb terms. This issue de- 
serves future attention. Very likely, our finding of large 
non-local Coulomb interaction Uqi in graphene general- 
izes to other two-dimensional materials. In narrow impu- 
rity bands or edge states of graphene, the Coulomb inter- 
action might in any case present the dominating energy 
scale and, thus, trigger many body instabilities including 
magnetism. 
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Here, e(k) is defined as the ratio of the bare and par- 
tially screened electron-electron interaction potentials 
with momentum transfer k, i.e., qi — > qi + k and 
q2 + k — > q2, averaged over the momenta qi and q2. 
The Wannier function index n is chosen to correspond to 
atoms in sublattice A. 
[37] The long wavelength behavior of e(k) determines the 
screening of the long range tails of the Coulomb inter- 
action. e(k) — > 1 for k — > corresponds to an unscreened 
1/r tail of the Coulomb interaction. 
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Matter 20, 494227 (2008). 
[39] Here, we consider a point charge in the middle of the 
film and evaluate the Coulomb potential in the middle 
of the film. Using the conventions and nomenclature of 
Ref. [3S| our situation corresponds to z = 0, d — h/2, 
k = 0, e 2 = e 3 = 1 and /3i 3 = P' 12 = (ei - l)/(ei + 1). 
Then, Eq. (4) of Ref. [3|| leads after division by the bare 
interaction (q/X) to e(k — A) as in our Eq. 



